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We solve the Hartree-Fock-Bogoliubov (HFB) equations for deformed, axially symmetric even- 
even nuclei in coordinate space on a 2-D lattice utilizing the Basis-Spline expansion method. Results 
are presented for the neutron-rich zirconium isotopes up to the two-neutron dripline. In particular, 
we calculate binding energies, two-neutron separation energies, normal densities and pairing den- 
sities, mean square radii, quadrupole moments, and pairing gaps. Very large prolate quadrupole 
deformations (fa = 0.42, 0.43, 0.47) are found for the 102 > 104 ' 112 Zr isotopes, in agreement with recent 
experimental data. We compare 2-D Basis-Spline lattice results with the results from a 2-D HFB 
code which uses a transformed harmonic oscillator basis. 

PACS numbers: 21.60.-n,21.60.Jz 
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I. INTRODUCTION 

In contrast to the well-understood behavior near the 
valley of stability, there are many open questions as we 
move towards the neutron dripline. In this exotic region 
of the nuclear chart, one expects to see several new phe- 
nomena 0: the neutron-matter distribution is very dif- 
fuse giving rise to the "neutron skins" and the "neutron 
halos" . One also expects new collective modes associ- 
ated with the neutron skin, e.g. the "scissors" vibrational 
mode. 

In current experiments, the neutron dripline has only 
been reached for relatively light nuclei. The second- 
generation radioactive ion beam (RIB) facilities currently 
under construction will allow us to explore nuclear struc- 
ture and astrophysics at the driplines of heavier nuclei. 

It is generally acknowledged that an accurate theoret- 
ical treatment of the pairing interaction is essential for 
a description of the exotic nuclei 0, 0- Besides large 
pairing correlations, the HFB calculations have to face 
the problem of an accurate description of the continuum 
states with a large spatial extent. All of these features 
represent major challenges for the numerical solution. 

There are various mean-field methods of solving the 
non-relativistic HFB equations. Generally, they can be 
divided into two categories: lattice methods and basis ex- 
pansion methods. In the lattice approach no region of the 
spatial lattice is favored over any other region: the well 
bound, weakly bound and (discretized) continuum states 
can be represented with the same accuracy. Currently, 
there are only two codes which solve the HFB problem on 
a lattice in coordinate space without any approximations 
using the quasiparticles: a 1-D (spherical) HFB code 3 
and the present 2-D (axially symmetric) code In 



the basis expansion method a wave function is expanded 
into the chosen basis functions: the harmonic oscillator 
basis (HO) [6j], the transformed harmonic oscillator basis 
(THO) 0, the HF orbitals HQ- There is also an alter- 
native method that is under development where the HFB 
equations are solved in the 3-D Cartesian mesh using the 
canonical-basis approach |lOj| . 

We solve the Hartree-Fock-Bogoliubov (HFB) equa- 
tions for deformed, axially symmetric nuclei in coordinate 
space on a 2-D lattice 0, Our computational tech- 
nique (the Basis-Spline collocation and Galerkin method) 
is particularly well suited to study the ground state prop- 
erties of nuclei near the driplines. It allows us to take into 
account high-energy continuum states up to an equiva- 
lent single-particle energy of 60 MeV or more. 

In this paper we study the ground state properties of 
neutron-rich zirconium nuclei (Z = 40) up to the two- 
neutron dripline. The isotope chain calculations start 
from the 102 Zr isotope up to the dripline nucleus, which 
turns out to be 122 Zr. For a long time, the A ~ 100 
region has been of interest to nuclear structure physi- 
cists as an area of competition between various coexisting 
nuclear shapes (well-deformed prolate, oblate, or spheri- 
cal) ^lj. The zirconium isotopes are known to possess a 
rapidly changing nuclear shape when the neutron num- 
ber changes from 56 to 60 [IJ- We find that a spheri- 
cal ground state shape is preferred over a prolate shape 
starting from the 114 Zr isotope up to the dripline nucleus 
122 Zr. In this paper, we compare Basis-Spline lattice re- 
sults with the HFB-2D-THO code [jj, which uses a trans- 
formed harmonic oscillator basis (THO), as well as with 
currently available experimental data. 
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II. SKYRME-HFB EQUATIONS IN 
COORDINATE SPACE 



III. NUMERICAL RESULTS AND 
COMPARISON WITH EXPERIMENTAL DATA 



A detailed description of the Basis-Spline lattice 
method has been published in references 0, Q ; we give 
here only a brief summary. In coordinate space repre- 
sentation, the HFB Hamiltonian and the quasiparticle 
wavefunctions depend on the distance vector r, spin pro- 
jection g = ±i, and isospin projection q = ±| (cor- 
responding to protons and neutrons, respectively). In 
the HFB formalism, there are two types of quasiparti- 
cle (bi-spinor) wavefunctions, (/)% and (f>2- In the present 
work, we utilize a Skyrme effective N-N interaction in the 
particle-hole channel while the p — p and h — h channels 
are described by a zero-range delta pairing force. For 
these types of effective interactions, the particle mean 
field Hamiltonian h and the pairing field Hamiltonian h 
are diagonal in isospin space and local in position space, 
and the HFB equations have the following structure 0: 



(h q - A) h« 
h q -(h q -X) 



= E a 



(1) 



where the Lagrange parameter A is introduced to yield 
the correct particle-number on average since the HFB 
wavefunction does not have a well defined particle- 
number. 

The quasiparticle energy spectrum is discrete for \E\ < 
—A and continuous for \E\ > —A For even-even nu- 
clei it is customary to solve the HFB equations for pos- 
itive quasiparticle energies and consider all negative en- 
ergy states as occupied in the HFB ground state. The 
quasiparticle wavefunctions determine the normal den- 
sity p q (r) and the pairing density p q (r) as follows 

p^) = E E <^(™) . (2) 

OO + 2 

-°» = E E ^J™) <t>l>°) ■ ( 3 ) 

E a >0a=-i 



Restricting ourselves to axially symmetric nuclei, we 
use cylindrical coordinates (r, z,(p). In the HFB lat- 
tice method, we introduce a 2-D grid (r Q , zp) with a — 
1, ...,N r and (3 = 1, ...,N Z . In radial direction, the grid 
spans the region from to r max . Because we want to 
be able to treat octupole shapes, we do not assume left- 
right symmetry in z-direction. Consequently, the grid 
extends from -z max to +z max . Typically, z max w r max 
and N z w 2 • N r . The HFB wavefunctions and operators 
are represented on the 2-D lattice by Basis-Spline expan- 
sion techniques. In practice, B-Splines of order M=9 are 
being used. Details about the B-Spline collocation and 
Galerkin methods are given in the Appendix. 



This section describes the numerical results for the 
even-even members of the zirconium (Z — 40) isotope 
chain up to the two-neutron dripline. Two theoretical 
approaches are presented and compared: a) a 2-D lattice 
method using B-Spline technology (hereafter referred to 
as HFB-2D-L ATTICE) , and b) an expansion in a 2-D 
harmonic oscillator (HFB-2D-HO) and transformed har- 
monic oscillator basis (HFB-2D-THO). We also compare 
theoretical results with experimental data whenever pos- 
sible. 

In all calculations we utilize the SLy4 Skyrme force 
|13| and a zero-range pairing force with a strength pa- 
rameter V = -187.1305 (HFB-2D-LATTICE code) and 
an equivalent single particle energy cutoff parameter of 
£max = 60 MeV. The pairing strength Vo has been ad- 
justed in both codes to reproduce the measured average 
neutron paring gap of 1.245 MeV in 120 Sn Q, as can be 
seen from Table [I] 



TABLE I: Results for foSn. Comparison of results ob- 
tained from the HFB-2D-LATTICE code (first row) with 
HFB-2D-THO results (second row) and the experiment (third 
row). The columns display binding energies (BE), intrinsic 
quadrupole moments for neutrons and protons (Qn,Q P ), rms- 
radii (r n ,r p ), average pairing gaps (A„,A P ), pairing energy 
for neutrons (PE n ), and Fermi levels (A ra , X p ). 



BE (MeV) 


Q„(fm 2 ) 


Q p (fm 2 ) 


r„ (fm) 


r p (fm) 


-1019.26 


0.29 


0.12 


4.725 


4.590 


-1018.22 


0.00 


0.00 


4.728 


4.593 


-1020.54 










A„ (MeV) 


A p (MeV) 


PE„(MeV) 


A ?l (MeV) 


A p (MeV) 


1.244999 


0.0 


-10.24 


-7.98 


-8.16 


1.245469 


0.0 


-10.26 


-7.99 


-11.13 


1.245 











In Table d we compare the results from the two HFB 
codes for the ^gSn isotope. All observables agree very 
well. The apparent "disagreement" in the proton Fermi 
level A p is really an artifact: the pairing gap vanishes 
at the magic proton number Z = 50 resulting in an 
ill-defined Fermi energy. The two codes use different 
prescriptions for calculating X p in the trivial case of no 
pairing. HFB+THO accepts the last occupied equiva- 
lent single-particle energy as A p in this no pairing case, 
whereas the Basis-Spline code takes the average of the 
last occupied and first unoccupied equivalent single- 
particle energy levels. 



A. Deformations, dripline and pairing properties 

Recently, triple-gamma coincidence experiments have 
been carried out with Gammasphere at LBNL 14] which 
have determined half-lives and quadrupole deformations 
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of the neutron-rich 102 < 104 Zr isotopes. The isotopes from 
that region are produced in the process of fission of 
transuranic elements and have been studied via 7-ray 
spectroscopy techniques. These medium-mass nuclei are 
among the most neutron-rich isotopes (N/Z ~ 1.6) for 
which spectroscopic data are available. Very large pro- 
late quadrupole deformations [fa = 0.43, 0.45) are found 
for the 102 ' 104 Zr isotopes. Furthermore, the laser spec- 
troscopy measurements [llj for the zirconium isotopes 
have yielded precise rms-radii in this region. Recently, 
an experiment has been carried out to measure the mass 
of the 104 Zr isotope [l^. It is therefore of a great in- 
terest to compare these data with the predictions of the 
self-consistent HFB mean field theory. 

In radial (r) direction, the lattice extends from 
0—15 fm, and in symmetry axis (z) direction from 
— 15, ...,+15 fm, with a lattice spacing of about 0.8 fm 
in the central region. Angular momentum projections 
ft = 1/2,3/2, ...,21/2 were taken into account. Figure^] 
shows the calculated two-neutron separation energies for 
the zirconium isotope chain. The two-neutron separation 
energy is defined as 

S 2n (Z,N) = E bmd (Z,N) - E bi nd(Z,N-2). (4) 

Note that in using this equation, all binding energies 
must be entered with a positive sign. The position 
of the two-neutron dripline is defined by the condition 
S2n(Z,N) = 0, and nuclei with negative two-neutron 
separation energy are unstable against the emission of 
two neutrons. As one can see, both codes (HFB-2D-THO 
and HFB-2D-LATTICE) predict the dripline nucleus to 
be at mass number A = 122. In addition, we observe a 
very good agreement between the two codes in the whole 
mass region of the isotope chain. We also give a com- 
parison with the latest available experimental data up to 
the isotope 110 Zr H^j. 
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FIG. 1: (Color online) Two-neutron separation energies for 
the neutron-rich zirconium isotopes. The dripline is located 
where the separation energy becomes zero. The 122 Zr isotope 
is the last stable nucleus against two neutron emission. 

In Fig. [3 we compare the intrinsic proton and neutron 
quadrupole moments calculated with the lattice code and 
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FIG. 2: (Color online) Intrinsic quadrupole moments for pro- 
tons and neutrons. 



the THO code. Available experimental data [14| are also 
given. Generally, we observe a nearly perfect agreement 
between the two codes as well as with the experiment. 
The deformations (for neutrons ) in terms of the defor- 
mation parameter fa for those nuclei, namely for the 
i02-ii2 Zr i so topes range from fa=0A2 to fa = 0.47. 
Both the Basis-Spline lattice code and the HFB-2D-THO 
code predict the 112 Zr isotope to have the largest ground 
state deformation. For mass numbers larger than 112, 
we observe spherical ground state shapes. Experimental 
deformations for protons are available for two isotopes, 
102 Zr and 104 Zr |l4|. Calculations agree with the experi- 
ment reasonably well and give fa values of 0.42,0.43 while 
the experiment predicts /3| 02 =0.42, /3| 04 =0.45. 
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FIG. 3: (Color online) Root-mean-square radii for the chain 
of zirconium isotopes. 

In Fig.Olwe compare the root-mean-square radii of pro- 
tons and neutrons predicted by the LATTICE code and 
the THO code. Both codes give nearly identical results 
for the whole isotope chain. Only one experimental data 
point is available, the proton rms radius of 102 Zr [llj . 
The experiment yields a proton rms radius of 4.54 fm 
while the HFB codes predict a value of 4.45 fm (HFB-2D- 



4 



LATTICE) and 4.46 fm (HFB-2D-TH0). The difference 
between theory and experiment is quite small, of order 
2%. We can clearly observe the presence of the neutron- 
skin manifested by the large differences between the neu- 
tron and proton rms radii for all of the isotopes in the 
chain. As expected the neutron-skin becomes "thicker" 
as we approach the dripline. Starting at the mass num- 
ber A=114 up to the dripline the nuclei prefer a spheri- 
cal ground state shape (Fig. |2J) which results in a sudden 
shrinking of the rms radius at A=114. 
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FIG. 4: (Color online) Average neutron pairing gap for the 
chain of zirconium isotopes. 
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FIG. 6: (Color online) Contour plots of the 110 Zr normal 
densities, for protons (left) and neutrons (right). Densities 
are shown as a function of the cylindrical coordinates (r, z), 
where z is the symmetry axis. 



B. Density studies 
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FIG. 5: (Color online) Average proton pairing gap for the 
chain of zirconium isotopes. 

Fig. 0| and Fig. depict the average pairing gaps for 
neutrons and protons. Generally, both HFB codes show 
the same trend for the pairing gaps as a function of neu- 
tron number; the agreement is noticeably better for neu- 
trons. The two HFB codes predict a small value of the 
neutron pairing gap for the 112 Zr isotope which on the 
other hand has the largest prolate deformation (Fig. |2J) 
among the calculated nuclei. Coincidentally, the dripline 
turns out to be at the neutron magic number (N=82) 
and, as expected, both codes yield a pairing gap of zero 
for the 122 Zr isotope. 



In this section we focus on the normal and pairing den- 
sities for the selected isotopes. In Fig. we show a con- 
tour plot of normal densities for protons and neutrons for 
the 110 Zr isotope. It is the last deformed isotope with a 
significant value of the pairing gap for neutrons, (Fig.^J), 
therefore it is possible to show plots of both normal and 
pairing densities. The results obtained for the neutron 
normal and pairing densities (Figs. and 0) clearly ex- 
hibit a large prolate deformation. The normal density 
for neutrons (Fig. is concentrated in the region that 
extends from fm to 2 fm in r-direction and from —5 fm 
to +5 fm in z-direction. Within this region, we find an 
enhancement in the neutron density with a shape that 
resembles the figure "eight". In comparing the neutron 
to the proton density, one notices that both the center of 
the nucleus and the surface is dominated by neutrons. 

The pairing density for neutrons in Fig. [7| shows a 
richer structure than the normal density. This quantity 
describes the probability of correlated nucleon pair for- 
mation with opposite spin projection, and it determines 
the pair transfer formfactor. We can see that most corre- 
lated pair formation takes place in the four closed shaped 
structured areas near the z-axis. We may conclude that 
neutrons dominate the pairing properties in the interior 
of this nucleus. The same applies to normal densities 
(Fig. 0), yet the difference between neutrons and pro- 
tons is more pronounced in case of the pairing densities. 
A graph depicting the single-particle energy spectrum of 
the pairing density for the 104 Zr isotope has been pub- 
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lished in Ref. 0. 

In Figs. 00 and [5] we show plots of normal densities as a 
function of the distance from the center, r = p 2 + z 2 . 
For a given value of r, the density is single- valued for a 
spherical nucleus and multi- valued for a deformed density 
distribution because in the latter case different combina- 
tions of lattice points Zi and pj give rise to the same 
r- value. 
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FIG. 7: (Color online) Contour plots of the 110 Zr pairing 
densities, for protons (left) and neutrons (right). The den- 
sities are shown as a function of the cylindrical coordinates 
(r, z), where z is the symmetry axis. 

In Fig.|Hlwe compare three different calculations of the 
neutron normal density for the most deformed 112 Zr iso- 
tope. The plot on a logarithmic scale shows that the 
density distribution predicted by the HFB-2D-THO and 
HFB-2D-LATTICE codes is deformed for almost all val- 
ues of the distance from the nuclear center, r. At very 
large distances the densities become less deformed since 
nuclear potentials go to zero and HFB equations lead to a 
spherical asymptotic solution. Fig. [8] also shows for com- 
parison the HFB-2D-HO result as an illustration of the 
shortcomings of the pure harmonic oscillator basis calcu- 
lations to reproduce density distributions asymptoticly 
at very large distances. One can see its too rapid de- 
cay beyond distances of about 12 fm. Clearly, the pure 
harmonic oscillator basis calculations cannot represent 
properly density asymptotic for nuclei close to the neu- 
tron drip line. 

Neutron and proton normal densities for the drip-linc 
nucleus 122 Zr are shown in Fig. [5] From the single- valued 
plot as a function of r = y/ p 2 + z 2 one can immediately 
conclude that both neutron and proton normal densities 
are spherical. Another striking feature is a strong neu- 
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FIG. 8: (Color online) Logarithmic plot of the normal neu- 
tron density for the most deformed isotope 112 Zr as a function 
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FIG. 9: (Color online) Linear plot of the normal neutron and 
proton density for the dripline nucleus 122 Zr as a function of 
distance r = y^p 2 + z 2 . Comparison between the HFB-2D- 
THO code and the HFB-2D-LATTICE code. 



tron enhancement at the center and a corresponding de- 
pletion in the proton density. Note also that the neutron 
density is substantially larger than the proton density for 
all values of r. 



IV. CONCLUSIONS 

In this paper, we have performed Skyrme-HFB calcu- 
lations in coordinate space for the neutron-rich zirconium 
isotopes up to the two-neutron dripline. We have calcu- 
lated the ground state properties (even-even nuclei) for 
the zirconium isotopes (Z=40), starting from N=62 up to 
the two-neutron dripline, which our HFB codes predict 
to be N=82. In particular, we have calculated the two- 
neutron separation energies, quadrupole moments, rms 
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radii, average pairing gaps and densities. In comparing 
HFB-2D-LATTICE theoretical calculations for the two- 
neutron separation energies (Fig. ^) with the HFB-2D- 
THO code we find the same results. Particularly, both 
codes predict the 122 Zr isotope to be the dripline nu- 
cleus. We find very large prolate deformations for the 
102-U22 r isotopes and a spherical ground state shape for 
the 114 - 122 Zr nuclei. The fa value for the most deformed 
nucleus 112 Zr in the calculated chain reaches an impres- 
sive value of 0.47. The root-mean-square radii clearly 
show the existence of a "neutron skin" in the neutron- 
rich zirconium isotopes. We can also observe a sudden 
shrinking of the rms-radius at A=114 due to a change of 
the prolate ground state deformation into the spherical 
shape. In section III.B. we studied normal and pairing 
densities. In particular, Figures |H1 and show the domi- 
nant nature of neutrons, where one can observe that both 
normal and pairing densities are dominated by neutrons. 
The density studies in Figure predict the presence of 
a sizable neutron skin in the spherical dripline nucleus 
122 Zr. 



APPENDIX: BASIS-SPLINE REPRESENTATION 
OF WAVEFUNCTIONS AND OPERATORS 



We consider an arbitrary (differential) operator equa- 
tion 

Of(x)-g(x)=0. (A.l) 

Special cases include eigenvalue equations of the 
HF/HFB type where O = h and g(x) = Ef(x). We 
assume that both f(x) and g(x) are well approximated 
by spline functions 

/»«/(*) = ]>>f(*K, 

z=l 
H 

g(x)*g(x) = ]Ti? 4 M (*)&*. (A.2) 
j=i 

Because the functions f(x) and g{x) are approximations 
to the exact functions fix) and g(x), the operator equa- 
tion will in general only be approximately fulfilled 

Of(x) - g(x) = R(x) . (A.3) 

The quantity R(x) is called the residual] it is a measure 
of the accuracy of the lattice representation. 



In this appendix, we discuss the representation of dif- 
ferential operators and wavefunctions in terms of Basis- 
Splines which is used in the HFB-2D-LATTICE code. 
There exists an extensive literature on Basis-Spline the- 
ory developed by mathematicians [23 ■ We have adapted 
some of these methods for the numerical solution of prob- 
lems in atomic and nuclear physics on a lattice. Ref. [2^ 
discusses the B-Spline collocation method, periodic and 
fixed boundary conditions, and the solution of various 1- 
D radial problems (Schrodinger, Poisson and Helmholtz 
equations) and the solution of 3-D Cartesian problems 
(Poisson equation). In a later paper, Umar et al. [24[ 
solved the HF and TDHF equations on a 3-D lattice. In 
1995, Wells et al. [25[ applied the B-Spline collocation 
method to the static and time-dependent Dirac equa- 
tion which eliminated the "fermion doubling problem" 
that one encounters with the finite-difference method. 
In 1996, Kegley et al. [2(| studied 2-D axially symmetric 
systems in cylindrical coordinates with the collocation 
method and also introduced the Basis-Spline Galerkin 
method. In the following, we compare both of these 
methods. 



1. HFB lattice representation 

For the lattice representation of wavefunctions and op- 
erators, we use a Basis-Spline method. Basis-Spline func- 
tions Bf 1 (x) are piecewise-continuous polynomials of or- 
der (M — 1). They represent generalizations of finite 
elements which are B-splines with M = 2. 



2. Basis-Spline collocation method 

In the collocation method, one minimizes the residual 
locally, i.e. one introduces a set of collocation (data) 
points x a and requires that the residual vanish exactly 
at these points 

R{x a ) = . (A.4) 

We multiply Ea. (|A.3|) from the left with S(x — x a ) and 
integrate over x, including a volume element weight func- 
tion v(x) in the integrals to emphasize that the formalism 
applies to arbitrary curvilinear coordinates. Most cases 
of interest are covered by a function of the form 

{0 Cartesian 
1 Polar . (A.5) 
2 Spherical 

Using the collocation condition Ea. (|A.4l) one obtains 

[Of(x)] Xa -g(x a ) = . (A.6) 

Inserting the B-Spline expansion (|A.2|I of the function 
f(x) this results in 

Y}OB] al a l -g{x a ) = , (A.7) 

i 

where we have introduced the shorthand notation 

[OB] m = {OBi{x)] Xa . (A.8) 
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We eliminate the expansion coefficients a 1 in Ea. (|A.7|) 
by introducing the function values at the lattice support 
points x a including both interior and boundary points 



f a = f(Xg) = ^ Bi{x a )a l = ^2 B a%a l 



(A.9) 



g a = g(x a ) = Y, B ^)b l = J2 B ^ bl ■ ( A - 10 ) 

i i 

The matrix B a i is, in general, rectangular. However, 
it can be made into a square matrix by adding cither 
periodic or fixed boundary conditions; this is described in 
detail in references |2a l25l l26| . In what follows, the new 
(square) B-Spline matrix with the boundary conditions 
added is also denoted by B a i, for simplicity of notation. 
Because the matrix B is now a square matrix it can be 
inverted to eliminate the expansion coefficients a 1 , b l 

a i = Y j B ia f a , b l = J2B* a g a , (A.ll) 

a a 

resulting in 

i 

Defining the collocation lattice representation of the op- 
erator O via 

oi = Y}° B ^ Blp i ( A - 13 ) 



we arrive at the desired lattice representation 

J2°^ = 9c- (A.14) 





3. Basis-Spline Galerkin method 

To derive the Galerkin representation, we multiply 
Ea. HA.3fl from the left with the spline function B k (x) 
and integrate over x 



v{x)dxB k {x)0 f (x) — J v(x)dxB k (x)g(x) — 
v(x)dxB k (x)R{x) . (A.15) 



Various schemes exist to minimize the residual function 
R(x); in the Galerkin method one requires that there be 
no overlap between the residual and an arbitrary B-splinc 
function 



v{x)dxB k {x)R(x) = . 



(A.16) 



This so called Galerkin condition amounts to a global re- 
duction of the residual. We apply the Galerkin condition 



to Ea. i|A.15[) and insert the B-Spline expansions for f(x) 
and g(x), Ea. HA.2[> . resulting in 



E 



E 



v{x)dxB k {x)OB l {x) 

v{x)dxB k (x)Bi (x) 

Defining the matrix elements 

O kl = J v{x)dxB k {x)OB l {x) , 
G ki = / v(x)dxB k (x)B l (x) . 



V = 0. (A.17) 



(A.18) 
(A.19) 



transforms the (differential) operator equation into a ma- 
trix x vector equation 



(A.20) 



which can be implemented on modern vector or paral- 
lel computers with high efficiency. The matrix G k i is 
sometimes referred to as the Gram matrix; it represents 
the non vanishing overlap integrals between different B- 
Splinc functions. This matrix possesses several highly 
desirable numerical properties; it is symmetric, banded, 
positive definite, and invertible for any reasonable place- 
ment of the knots and any spline order 



(A.21) 



We use again the relations l|A.ll() to eliminate the ex- 
pansion coefficients a 1 and b % in Ea. l|A.20|) which results 
in the matrix equation on the Galerkin lattice 







(A.22) 



with the differential operator definition on the Galerkin 
lattice 



O^J2 B ^ Gll °^B kp 

ijk 



(A.23) 



We are also extending our previous B-spline work to 
include nonlinear grids. Use of a nonlinear lattice should 
be most useful for loosely bound systems near the proton 
or neutron drip lines. Non-Cartesian coordinates necessi- 
tate the use of fixed endpoint boundary conditions; much 
effort has been directed toward improving the treatment 
of these boundaries. 



4. 2-D lattice representation of HFB wavefunctions 
and Hamiltonian 



For a given angular momentum projection quantum 
number f2, we solve the eigenvalue problem on a 2-D 
grid (r a , zp) where a = 1, m and f3 = 1, n. 
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The four components of the spinor wavefunction 
■0( 4 )(r, z) are represented on the 2-D lattice by a product 
of Basis Spline functions Bi{x) evaluated at the lattice 
support points. For example, the spin-up component of 
the wavefunction (j) 2 is represented in the form 

tfirczp, T) = Yl <M Bf{zp) U? . (A.24) 
hi 

Hence, the four-dimensional spinor wavefunction in coor- 
dinate space ip( 4 \r, z) becomes an array ip(N) of length 
N = 4 • m ■ n. 

For the lattice representation of the HFB Hamiltonian, 
we use a hybrid method in which derivative operators are 
constructed using the Galerkin method; this amounts to 
a global error reduction. Local potentials are represented 
by the B-Spline collocation method (local error reduc- 
tion). 

The lattice representation transforms the differential 
operator equation into a matrix x vector equation 



N 



(A.25) 



r - r 



1 



|r| r 



r v- 



1 ( , -\ 2 1 
2l r ' V ) r 

4 1 

- + ..., 

r 



with the following boundary conditions 



<t>c{p 



%max) 



o, 

. 



(A.27) 



The differential operators are represented on the lattice 
using the B-Spline Galerkin method. 
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5. 2-D B-Spline Poisson Solver 

In the current version of our HFB-2D-LATTICE code 
we have solved the Poisson equation in cylindrical coor- 
dinates (p, z) 



d 2 ^(p,z) , cj>{p,z) , d 2 4>{p,z) 



dp 2 



by applying a large distance expansion of the Coulomb 
potential (|r| > |r'|) 



; (r) = e 



/ 



(A.26) 
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